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ABSTRACT 

It has been realised only recently that TeV emission from blazars can significantly 
heat the intergalactic medium (IGM) by pair-producing high-energy electrons and 
positrons, which in turn excite vigorous plasma instabilities, leading to a local dissi- 
pation of the pairs' kinetic energy. In this work, we use cosmological hydro dynamical 
simulations to model the impact of this blazar heating on the Lyman-a forest at in- 
termediate redshifts {z ~ 2 — 3). We find that blazar heating produces an inverted 
temperature-density relation in the IGM and naturally resolves many of the problems 
present in previous simulations of the forest that included photoionisation heating 
alone. In particular, our simulations with blazar heating simultaneously reproduce the 
observed effective optical depth and temperature as a function of redshift, the observed 
probability distribution functions (PDFs) of the transmitted fiux, and the observed 
fiux power spectra, over the full redshift range 2 < 2: < 3 analysed here. Additionally, 
by deblending the absorption features of Lyman-a spectra into a sum of thermally 
broadened individual lines, we find superb agreement with the observed lower cutoff 
of the line- width distribution and abundances of neutral hydrogen column densities 
per unit redshift. Using the most recent constraints on the cosmic ultraviolet (UV) 
background, this excellent agreement with observations does not require rescaling the 
amplitude of the UV background; a procedure that was routinely used in the past 
to match the observed level of transmitted fiux. We also show that our blazar-heated 
model matches the data better than standard simulations even when such a rescaling 
is allowed. This concordance between Lyman-a data and simulation results, which 
are based on the most recent cosmological parameters, also suggests that the inclu- 
sion of blazar heating alleviates previous tensions on constraints for a^ derived from 
Lyman-a measurements and other cosmological data. Finally, we show that blazar 
heating dramatically alters the volume- weighted temperature PDF, implying an im- 
portant change of the strengths of structure formation shocks (and thereby possibly 
particle acceleration in these shocks). The density PDF is also modified, suggesting 
that blazar heating may have interesting effects on structure formation, particularly 
on the smallest galaxies. 

Key words: cosmology: theory - methods: numerical — intergalactic medium 



1 INTRODUCTION 

The thermal history of the intergalactic medium (IGM) is 
of key importance for understanding cosmological structure 
formation since it determines the initial thermodynamic con- 
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ditions for the formation of collapsed objects ranging from 
dwarf galaxies up to the scale of galaxy clusters. Simultane- 
ously, the IGM acts as a cosmic heat reservoir that records 
reionisation processes in the Universe that inject substantial 
amounts of energy into the IGM on comparably short cos- 
mological time scales. Detailed measurements of the thermal 
history of the IGM may therefore reveal the main contrib- 



2 E. Puchwein, C. Pfrommer, V. Springel, A. E. Broderick & P. Chang 



utors to the UV background, e.g., cosmic star formation or 
quasars, which are thought to be primarily responsible for 
these reionisation events. 

Until recently, it was generally thought that all the main 
sources for heating the IGM had been identified and only the 
detailed interplay of these processes remained t o be under- 



stood. However, in a recent series of 



^_^__^^,. - — , papers ( Broderick et al.l 

I2OIII : IChang et all I2OIII : IPfrommer etaP l201lh . a novel 
heating process was proposed that taps efficiently into the 
energy reservoir of TeV blazars which are ultimately pow- 
ered by accretion onto their central super-massive black 
holes (or more generally the engines of active galactic nu- 
clei, AGN). A subset of blazars emits the bulk of their energy 
in TeV 7-ray radiation {E > 100 GeV) to which the Uni- 
verse is opaque; i.e., these energetic gamma rays necessar- 
ily annihilate upon the photons comprising the extragalac- 
tic background light (EBL), producing an ultrarelativistic 
population of electron-positron pairs. Typical gamma-ray 
mean free paths range from 30 Mpc to 1 Gpc, depending 
upon 7-ray energy and source redshift. The pairs produced 
by TeV 7-ray s have typical Lorentz factors of 10^ — 10^. 
Despite its dilute nature, a beam of these ultra-relativistic 
pairs is susceptible to plasma beam instabilities while prop- 
agating through the IGM. In particular, the "oblique" insta- 
bility, a more virulent and robust cousin of the commonly 
discussed Weibel and two-stream instabilities, appears to 
grow o ut of the pair beam anisotropy as a source of free 
energv JBret et al.ll2004l. [20051 : lBretil2009l : iBret et aLllJoTol : 
iLemoine fc Pellet ierll20 id ). Though the nonlinear evolution 
of these instabilities is unknown, they may saturate at a level 
such that they dissipate the kinetic energy of the beam at a 
rate comparable to their linear growth rate, thereby heating 
the IG M locally, following the suggestion of iBroderick et all 

dmll). 

If this scenario or an analogous, similarly efficient mech- 
anism operates in practise, it necessarily suppresses the in- 
verse Compton cooling of the pairs, which acts on a longer 
timescale. This naturally explains the absence of an inverse 
Compton bump in observed TeV-blaz ar spectra without in- 
voking an intergalactic magnetic field ([Broderick et al.ll201U 
and references therein). At the same time, it allows for a 
redshift evolution of TeV blazars that is identical to that of 
quasars without overproducing the extragalactic 7-ray back- 
ground (EGRB). In fact, for plausible parameters of TeV- 
blazar spectra, it is possible to expla in the high-energy part 
of the EGRB (IBroderick et al.ll201lh . This is in stark con- 
trast to the previous picture which assumed that the TeV 
7-ray emission is reprocessed into the GeV regime by means 
of inverse Compton processes. In this picture, in order to 
not o verproduce the EGRB observed by Fermi (|Abdo et al.l 
I2OIOI), the spectral redistribution of TeV 7-ray energy into 
the GeV regime strongly constrains the evolution of the lu- 
minosity density of TeV blazars. In these studies it was found 
that TeV blazars cannot exhibit the dramatic rise in number 
seen in the quasar distribution by z ~ 1-2, i.e., the comoving 
num ber of blazars must have remained essentially fixed (see , 
e.g.jNarumoto fc Tot ani' 2006 ^ iKneiske Sz MannheimI I2OO8I : 
llnoue fc Totanill2009l : .Vent era 1 20 id ): at odds with both the 
large-scale mass assembly history of the Universe as traced 
by the star formation history, and with the luminosity his- 
tory of similarly accreting systems, e.g., the quasar luminos- 
ity density. 



Integrating over the energy flux per mean free path 
of all known TeV blazars yields a luminosity density, or 
equivalent ly a local heating rate, that dominates that of 
photoheating by more than an order of magnitude at the 
present epoch at mea n density, after acc ounting for incom- 
pleteness corrections (^ Chang et al.ll201lh . As demonstrated 
in Broderick et al. (201l|), the local TeV blazar luminosity 
function is consisten t with a scaled vers ion of the quasar 
luminosity function ([Hopkins et al.l 120071 ). thus the conser- 
vative assumption is that they evolve similarly, presumably 
due to the same underlying accretion physics. With this as- 
sumption, blazar heating starts to become important around 
the epoch of He 11 reionisation (z ^ 3.5). By 2; = 0, the blazar 
heating rate is more than an order of magnitude larger in 
comparison with the photoheating rate alone. 

Because the Universe is opaque to the propagation of 
TeV photons over cosmological distances, the redshift evo- 
lution of the TeV blazar luminosity density is observation- 
ally inaccessible for z > 0.5. Also, the exact value of the 
blazar heating rate is subject to an uncertain incomplete- 
ness correction factor, which accounts for the incomplete 
sky coverage of current TeV instruments, the duty cycle of 
TeV blazars, their spectral variability, variations of the TeV- 
blazar redshift evolution, and additional TeV sources that 
may contribute to the plasma instability heating. Hence we 
can try to invert the problem and investigate how we can 
use the IGM as a cosmic calorimeter to probe the physics 
of TeV 7-rays and the global properties of TeV blazars by 
means of the energy that is deposited through plasma insta- 
bilities (or similarly efficient dissipation processes of the TeV 
emission from blazars). In this manner, the thermal history 
of the IGM provides an analogous argument to that of lSoltanI 
([l982l ). yielding an important constraint on the redshift evo - 
lution of the blazar luminosity density ([ Chang et al.ll201lh . 
Hence we turn to precision measurements o f the Lyman- 
a forest at z ^ 2 — 3 (e.g., IViel et al.l l2009l l in this work, 
with the goal to clarify the possible signatures of such an 
additional heating process, which is characterised by some 
unusual properties as we will discuss now. 

Since the number densities of EBL photons and of TeV 
blazars are nearly homogeneous on cosmological scales, so is 
the resulting pair density. Hence, the implied heating rate 
is also homogeneous, i.e., the volumetric blazar heating rate 
is expected to be spatially uniform and independent of IGM 
density. This implies that the energy deposited per baryon is 
substantially larger in more tenuous parts of the Universe, 
such that underdense regions experience a larger temper- 
ature increase, producing an inverted temperature-density 
(T-p) relation with an asymptotic scaling of T oc 1/p for 
densities much lower than the cosmic mean. Interestingly, 
detailed studies of high-resolution Lyman-a forest spectra 
led to the conjecture t hat such an inverted T-p relation i n 
fact exists at ^ = 2 - 3 (iBolton et al.ll2008l : IViel et aLllJOOol ) , 
something which has proven hard to a rrange within the con- 
text of standard reioni sation models ([McQuinn et al.ll2009l : 
JBolton fc Beckerll2009h . Related to that, when deblending 
the absorption features of Lyman-a spectra into a sum of 
thermally broadened individual lines and measuring their 
strengths and widths, simulations produce line widths that 
are significantly narrower than observed, implying an associ- 
ated inconsistency b etween the simulat e d and measured line 
width distribution ([Brvan et al.l Il999l : JBrvan fc Machacekl 
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|2QQQ| : iTheuns et alJiToOOl : iMachacek et alJlJOOOl ). This find- 
ing is largely independent of the amount of small-scale 
power, but does depend strongly on the IGM temperature; 
heuristically increasing the temperature of the IGM can 
henc e reproduce the obse rved line- width distribution rather 
weh (iTheunset al.ll200Ql ). 

In this work, we perform state-of-the-art simulations of 
the Lyman-a forest that include different variants of blazar 
heating (to address the uncertainty of the overall heating 
rate), and we compare the predictions to data in observa- 
tional space. In this way, we obtain important information 
on the validity of the blazar heating scenario and we can 
elucidate whether it helps to resolve current problems in 
understanding some aspects of the IGM observations. Af- 
ter describing our methods for simulating the Lyman-a for- 
est in Sect. [2l we present our results in Sect. [S] We focus 
on the thermal history and the T-p relation of the IGM in 
Sect. 13.11 the effective optical depths and photoionisation 
rates in Sect. [321 the transmitted flux probability distribu- 
tion functions and power spectra in Sect. 13.31 and the Voigt 
profile fitting of the Lyman-a forest in Sect. 13.41 Finally, 
we give an outlook on implications for cosmological param- 
eter estimation when including the physics presented in this 
work in Sect.fH and conclude in Sect.O 



2 METHODS 

2.1 The cosmological simulations 

Throughout this paper, we make use of a set of cosmolog- 
ical hydrodynamical simulations that were performed with 
the TreePM-SPH simulation code P-Gadget-3, which is 
an u pdated and extended version of Gadget- 2 (jSpringell 
I2OO5I). The simulations a dopt the WMAP 7-year cosmol- 
ogy (|Komatsu et al.ll2Qllh with Qm = 0.272, Qa = 0.728, 
r^B = 0.0465, h = 0.704, and ag = 0.809. Initial condi- 
tions with 2 X 128^ and 2 x 384^ particles (half of them gas, 
the other half dark matter) in a cosmological box with a co- 
moving side length of 15 /i~^ Mpc were evolved from redshift 
^^start = 100 to redshift Zend = or Zend = 1-9, respectively. 

Our analysis of the simulated Lyman-a forest is solely 
based on the higher resolution simulations with a fixed co- 
moving gravitational softening of e = 1.5 /i~^ kpc (Plummer- 
equivalent) and dark matter and gas particle masses of 
Mdm = 3.7 X IO^H-^Mq and Mgas = 7.5 x lO^h'^MQ 
respectively, while the lower resolution runs with a three 
times larger softening and 27 times larger particle masses 
were only used to study the overall thermal history of the 
universe G 

All of our runs follow hydrodynamics with an 
entropy- conserving formulation of smoothed-particle hydro- 
dynamics (SPH) fsee Springel fc Hernquist 2002 ), and ac- 
count for radiative cooling, star formation, and super- 
novae feedback with a sub- resolution model described in 
ISpringel fc HernquistI (|2003l V Photoheating was included 



under the assumption of ionisation equilibrium in the pres- 
ence o f an external UV background field as in iKatz et al.l 
(| 19961 ). however, using a state-of-the-art model for the 
time evolution of the cosm ic UV background (based on 
iFaucher-Giguere et al.ll2009l . hereafter FG'09, corresponding 
to their Table 2). This UV backgro und m odel ag rees well 
with observational constraints from iFaucher-Giguere et al.l 
(I2OO8I ) (hereafter FG'08). 

Typically, when photoheating is followed in such a 
framework, hydrogen (H) and helium i (He i) reionisation 
happen very rapidly, which is problematic for two reasons: 
first, the assumption of ionisation equilibrium may not be 
satisfied during the reionisation transition itself, and sec- 
ond, the time steps in the simulation (which are limited by 
the Courant condition) may be too large to yield the correct 
overall temperature increase even if the equilibrium assump- 
tion were still valid, as the sharp peak in the heating rate 
may not be sampled well. This usually results in a too in- 
efficient heating. We circumvent both problems by injecting 
the overall heat input expected for a sudden H and He i 
reionisation by hand. For this, we need to assume values for 
the average excess energy per H and He I ionisation event. 
We chose values that are consistent with the FG'09 model 
at ^ = 10, namely, Eu = 5.8 eV and ^Hei = 9.6 eV. This 
results in a temperature increase of ATH,Hei = 2.35 x 10^ K. 

One possible way in which to apply such a temperature 
boost would be to turn off the standard photoheating that 
assumes ionisation equilibrium with the UV background un- 
til z :^ 10. We could then apply the temperature boost, 
ATH,Hei, to the SPH particles and subsequently turn on the 
standard photoheating to follow He 11 reionisation. At that 
time, our assumed UV background would be strong enough 
to keep H and He i ionised. 

However, we adopted a different approach that allows us 
to bracket the uncertainties involved in the treatment of H 
and He i reionisation. This involved performing two simula- 
tion runs where we did not disable the standard photoheat- 
ing. In one of them we applied the temperature boost shortly 
after the reionisation by the UV background. This procedure 
somewhat overpredicts the heating since the standard pho- 
toheating implementation, despite being inefficient, will also 
have injected some energy. In the other run, we injected the 
energy slightly before the UV background reionises H and 
He I. In this case, even after the injection, a fraction of the H 
and He will still be considered neutral given the assumption 
of ionisation equilibrium with the UV background. Thus, 
atomic cooling will remove part of the injected energy from 
the IGM until the UV flux has increased sufficiently to fully 
ionise H and He I. Hence, in the second approach, we un- 
derestimate the temperature after H and He i reionisation. 
Together, the two runs bracket the actual heat input. The 
thermal history of the IGM in these runs is discussed in more 
detail in Sect. 13. Il 



^ Results from the 2 x 128"^ runs are only shown in Figs. [2] and 
3] In the 2 x 384^ simula tions, the Lyman -o; forest properties are 
well converged (see e.g., IViel et al.|[2004l , for a discussion of the 
convergence of power spectra of the transmitted flux). 



Finally, we would like to point out that He 11 reioni- 
sation happens much more gradually in the FG'09 model, 
with the bulk of the He 11 being reionised between redshifts 
5 and 3. Hence, we are confident that the standard photo- 
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heating algorithm accurately follows the He ii reionisation 
and yields the correct heating rates for the FG'09 modellj 



2.2 Blazar heating 

The main focus of this work is an investigation of the impact 
of blazar heating on the properties of the IGM and its obser- 
vational signatures, particularly the Lyman- a forest. Thus, 
we performed four simulations at both resolutions, one with- 
out blazar heating and three with blazar heating at different 
levels of efficiency. All four runs were started from identical 
initial conditions, so they are directly comparable, allowing 
us to isolate and cleanly study the effects of blazar heating. 
In the simulations, we turn on blazar heating at redshift 
z = 5 and adopt the redshift evolution of the b lazar heating 
rate p er co-moving volume Q{z) suggested bv I Chang et al.l 
(l201lh for ^ < 5, i. e. 



simulation model 



Q{z = 0) 
[eV Gyr"-*^ cm 



^sys 



■0.512 X [(1 + ; 



l]+2.27x[(l + z)-l], (1) 



where our assumed volumetric heating rates at the present 
epoch, Q{z — 0), are summarised in Table [T] Our weak 
blazar heating model is sligh t ly les s efficient than the stan- 
dard model in IChang et al.l (|201ll ). while their optimistic 
model falls in between our intermediate and strong heat- 
ing. For completeness, we also list the correction factors ?7sys 
th at correspond to o ur heating rates. ?7sys (defined in Eq. 16 
of lChang et al.ll201lh is a parameter of order unity that cor- 
rects for remaining systematic uncertainties in the predic- 
tion of the heating rate due to the spectral variability of 
TeV blazars, their uncertain redshift evolution, and possible 
additional TeV sources that may contribute to the plasma 
instability heating. 

Based on the rate given by Eq. ([1]), the required heat 
injection for each SPH particle is implemented in the simu- 
lation code by accordingly updating the internal energy of 
all active SPH-particles at each time- step. The required esti- 
mate of the volume, V = mspu/pspu, occupied by a particle 
is computed as ratio of particle mass, mspH, to SPH density 
estimate, psPHU 



2.3 Simulating the Lyman-a forest 

In order to be able to compare the IGM predicted by our 
simulations to observations of the Lyman-a forest, we need 
to compute synthetic Lyman-a forest absorption spectra. 
We do this by selecting a large number of randomly placed 
lines of sight through individual outputs of the simulation 
box, along directions parallel to one of the coordinate axes 
(randomly selected among x, ^, and z), and each represented 



^ There is, however, still some uncertainty in the He ll heating 
rate due to assuming optical-thin photoheating. 
^ In SPH simulations, the particle volumes obtained in this way 
typically do not add up exactly to the total volume of the simu- 
lation box. The sum of the volumes is about 5% low in our runs, 
which potentially biases the heating rate low by the same frac- 
tion. We do not attempt to correct for this as this bias is small 
compared to the uncertainty in the blazar-heating rate. 



no blazar heating 
weak blazar heating 
intermediate blazar heating 
strong blazar heating 




5.8 X 10-s 
1.08 X 10-'^ 
1.62 X 10-'^ 




0.66 
1.23 

1.85 



Table 1. Assumed volumetric blazar heating rates at z = in our 
different simulations. The redshift evolution is in all cases given 
by Eq. ([T} . Also listed are the values of ?7sys corresponding to our 
heating rates. 



by 2048 pixels. Next, we compute the neutral hydrogen den- 
sity, temperature and velocity fields of the IGM along these 
lines of sight by adding up the density contributions and av- 
eraging the temperatures and velocities of all SPH particles 
whose smoothing lengths are intersected. The calculation 
of the spectra then accounts for Doppler-shifts due to bulk 
flows of the gas along the line of sight as well as for ther- 
mal broadening of the Lym an- a resonance in the rest frame 
of absorbing gas (see e.g. iBolton fc Haehneltl l2007l V This 
procedure yields the optical depth, r, for Lyman-a absorp- 
tion as a function of velocity offset along each line of sight. 
This can then be easily converted into normalised transmit- 
ted flux, F = e~^, as a function of wavelength or redshift. 
Examples of simulated absorption spectra, obtained in this 
way, are shown in Fig. [U 



3 RESULTS 

3.1 The thermal history and the 

temperature-density relation of the 
intergalactic medium 

In most of the models for the thermal history of the 
IGM found in the literature presently, hydrogen and he- 
lium reionisation injects the majority of the energy. Recent 
observations of the cosmic microwave background suggest 
that hydrogen was r eionised around redshift Zreion ~ 10 
(jKomatsu et al.ll201lh O The UV background adopted in our 
simulations is consistent with this finding and reionises hy- 
drogen quickly at z = 10. 

As already discussed in Sect. 12. 1[ we have performed 
two runs with different treatments of H and He i reionisa- 
tion. They were designed to bracket the actual heat input 
during this epoch and were both performed without blazar 
heating. The thermal history of the IGM at the mean den- 
sity is illustrated for both of these simulations in Fig. [2l 
Shortly after H and He i reionisation, the temperatures are 
very different. However, adiabatic cooling due to the Hub- 
ble expansion, as well as by inverse- Compt on scattering off 
the CMB subsequently removes most of the injected en- 
ergy. By 2; = 5, this loss- of -memory effect has made the 
temperatures in both runs virtually identical, since we em- 
ploy a model with sufficiently early H and He i reionisation 



^ Their constraint ^^reion = 10-6 ±1.2 is, however, weakened if 
a gradual reionisation or multiple epochs of H reionisation are 
considered. 
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velocity [km s ^] 

Figure 1. The Lyman-o; forest in simulations with intermediate and without blazar heating at redshift z = 3. More precisely, the top panel 
shows the transmitted flux fraction as a function of velocity shift for three randomly selected lines-of-sight through the simulation box. 
The three lines-of-sight were stitched together at suitable positions (indicated by the vertical dotted lines) to allow a better comparison 
to real spectra. Clearly, Lyman-o; absorption lines are signiflcantly deeper in the run without blazar heating. The bottom panel shows 
the increase of the transmitted flux fraction due to blazars. The spectra were computed based on the FG'09 UV background. 



(jHui fc Haimanll2QQ3l ). Overall, this shows that the proper- 
ties of the IGM at z < 5 are comparatively insensitive to the 
treatment of H and He I reionisatiorjj, essentially because 
this happens so early, whereas heat input at later times, 
e.g., due to He ii reionisation, is more important in compar- 
ison. In the following, we only use simulations in which the 
ATH,Hei-boost was applied slightly before z = 10. 

For z < 5, blazar heating becomes the dominant heating 
mechanism around and below mean densities in all simula- 
tions that include this effect. As a consequence, the IGM 
temperatures rise significantly in these runs and reach val- 
ues that are about an order of magnitude larger than the 
temperature in the reference run without blazar heating. 
This dramatic effect is also illustr a ted in Fig. [2l 

As suggested in IChang et al.l (|201lh , blazar heating is 
implemented in our simulations using a volumetric heating 
rate that is independent of density. Consequently, less dense 
regions receive more heating energy per unit gas mass, caus- 
ing a larger temperature increase there. This results in an 
inverted T — p relation of the IGM, as can be clearly seen 
in Fig. [3l which shows the distribution of gas particles in 
the density and temperature plane at redshift z = 3, both 
for the simulation with intermediate blazar heating and for 
the reference run without blazar heating. In the latter case. 



^ Again, this statement may be weakened if a gradual reionisation 
or multiple epochs of H and He I reionisation are assumed. 



most of the low- density IGM follows a power- law T-p rela- 
tion, with a slope that is set by photoheating and adiabatic 
expansion and compression of the IGM due to structure for- 
mation. When including blazar heating, the T-p relation of 
the low-density IGM no longer has a constant slope. Instead, 
we find an inverted T-p relation in underdense regions, which 
turns over in an T-p relation with a positive slope at higher 
densities where blazar heating becomes progressively less ef- 
ficient. At 2; = 3, this turn-over occurs roughly at mean den- 
sity. It shifts, however, to somewhat larger overdensities at 
later times. 



IViel et ahl (|2009l l found that the transmitted flux dis- 
tributions of the observed Lyman-a forest are, indeed, best 
fit by an inverted T-p relation. Their best fits based on the 
ranges F — 0.1 — 0.8 and F — — 0.9 for the transmitted 
flux fraction are also indicted in the figure. The comparison 
is shown only for the range of overdensities that contributes 
most to these transmitted flux intervals. More importantly, 

I 1 r j x 

the inverted T-p relations inferred bv lViel et al.l (|2009l ) are 
clearly not found in the simulation without blazar heating, 
while the simulation with intermediate blazar heating is in 
agreement with these results. 



Also shown in Fig. [3] are contours of the neutral hydro- 
gen fraction, which where computed under the assumption 
of ionisation equilibrium with the UV background. Obvi- 
ously, the inverted T-p relation reaches lower neutral hydro- 
gen fractions for low densities. Hence, we expect significant 
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no blazar heating intermediate blazar heating 




3 -2 

logio(p/{P» 

Figure 3. The gas density-temperature distribution in simulations without blazar heating (left panel) and with intermediate blazar 
heating (right panel) at redshift z = 3. The gas mass per pixel Mpjx is colour-coded. The colour scale spans 5 orders of magnitude. 
Th e contours i ndica te the neutral hydrogen fractions log]^o(^Hl/^H) in the simulations. Also shown are the IGM T-p relations obtained 
by IViel et al.l (|2009l ) based on an analysis of the Lyman-a forest transmitted flux distribution in two transmitted flux fraction ranges, 
F = 0.1 — 0.8 and F = — 0.9, respectively, whose difference indicates the uncertainty in the data. 



o 



no blazar heating 
intermediate b. h. 
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Figure 4. Volume-weighted probability distribution functions of temperature (left panel) and density (right panel) of the IGM. Results 
are shown for redshifts z between and 5.3, and for simulations without blazar heating (dotted) and with intermediate blazar heating 
(solid). In the latter run, the temperature distribution is shifted towards larger values once blazar heating becomes efficient for z < 5. 
Without blazar heating, the temperature distribution evolves towards lower values once adiabatic cooling by the Hubble expansion 
exceeds heating by He ll-reionisation at z < 3.5. With blazar heating, Hubble cooling starts to dominate only at z < 0.5. In both runs, 
the density distribution shifts towards lower values (in units of the mean baryon density (p)) as a consequence of structure formation. It 
can also be clearly seen that blazar heating has a noticeable effect on the density distribution for z ^ 2. 
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Figure 2. IGM temperatures in simulations with and without 
blazar heating as a function of redshift z. Shown are the median 
temperatures of the gas particles with densities between 0.95 and 
1.05 times the mean cosmic baryon density. At low redshift, blazar 
heating strongly boosts IGM temperatures. For simulations with- 
out blazar heating, the temperature evolutions are shown for two 
different methods of treating hydrogen reionisation. We show that 
the IGM temperature at redshifts 2; < 5 does not depend on the 
implementation of hydrogen reionisation. The x-axis was chosen 
to be linear in \o^^^{\ + z). The no blazar heating curve is some- 
what difficult to identify since it is identical with the blazar heat- 
ing curves for z > 5 and almost identical with the no blazar 
heating, modified H-reionization curve for z < 5. 



changes in the Lyman-a forest. This will be explored in more 
detail in the following sections. 

Overall, Fig.[3]shows that the distribution of IGM mass 
in temperature and density is significantly altered by blazar 
heating. The changes are, however, even more dramatic in 
the distribution of the cosmic volume in temperature and 
density. This is illustrated in Fig. S] At redshift z = 5.3, 
before blazar heating is turned on, both simulations natu- 
rally have the same temperature and density distributions. 
Subsequently, due to ongoing structure formation, most of 
the volume evolves towards lower density compared to the 
mean density. Such underdense regions evolve, however, very 
differently in runs with and without blazar heating. In the 
latter case, they cool adiabatically by the Hubble expansion 
once the heating by He 11 reionisation becomes subdominant 
at 2; < 3.5. In contrast, blazars heat these regions efficiently 
in the former case. Together, this results in fundamentally 
different evolutions of the temperature distributions in runs 
with and without blazar heating. By redshift 2; = 0, the 
IGM that permeates most of the volume astonishingly is 
two orders of magnitude hotter when blazar heating is in- 
cluded. At low redshift, the differences in the density distri- 
bution also become noticeable. In particular, there is less vol- 
ume above the mean density, which likely has implications 
for lat e-forming structures as suggested bv lPfrommer et al.l 
(|2Qllh . 

In Fig. O we directly compare the temperature in our 
simulations to tempe rature measuremen ts based on Lyman- 
a forest observations (jBecker et al.ll201lh . The observational 



constraints were obtained by relating the curvature of ob- 
served Lyman-a forest spectra to the IGM temperature at a 
specific redshift-dependent overdensity A{z). A{z) increases 
from A '-^ 1 at z = 5 to ~ 6 at 2; = 2. We computed the IGM 
temperatures in our simulation at the same overdensities to 
make them directly comparable. 

However, even at a fixed overdensity there is some am- 
biguity in the temperature definition. This can also be seen 
in Fig. [3l which shows that the temperature distribution at 
a fixed overdensity between 3 and 6 is skewed and has a 
long tail towards high temperatures. Thus, the mean, the 
media n, and the mode of the temperature distribution will 
differ. iBecker et all (|2Qllh define the temperature T(A) by 
fitting a power law T-p relation near the mean density and 
evaluate this fit at overdensity A. In runs without blazar 
heating, the temperature obtained in this way closely tracks 
the mode of the mass- weighted temperature distribution at 
overdensity A, which is the temperature where the distribu- 
tion attains its maximum value. Thus, we use this temper- 
ature definition in Fig. [5] 

For 2; < 3, the simulation with intermediate blazar heat- 
ing is in excellent agreement with the data, while the run 
without blazar heating has too low temperatures T(A). 

For z > 3, where the temperature increase due to blazar 
heating becomes progressively unimportant, the IGM tem- 
peratures in all our simul ations, includi i ig the one without 
blazar heating, exceed the lBecker et al.l (|201lh values. This 
indicates tension between the early He II reionisation history 
implied by the FG'09 UV background model an d the Lyman - 
a forest temperature measurements by Becker et al.l (|201lh . 
We would like to point out, however, that our blazar heat- 
ing simulations are also in good agreement with temperature 
measurements based on the small-scale Lyman-a flux po wer 
spectrum (e.g., IZaldarriaga et aUl200ll : iLidz et aUl2009l l in 
the range 2 < z < 3, while our simulated temperatures are 
smaller than these measured values for z > 3.5. For z > 3.5, 
the simulated temperatures in all our runs, thus, fall in be- 
tween constraints based on the flux curvuture statistics and 
the small-scale power spectrum, while being formaly incon- 
sistent with both at the 1 — a level. This points to an uniden- 
tified systematics in either one or both of the methods for 
z > 3.5, where the absorbed fraction of the spectrum in- 
creases considerably, possibly suggesting an observational 
complication in measuring the temperature of the IGM at 
these redshifts. 

Overall, the thermal histories of our simulations with 
and without blazar heating are remarkably different. In the 
next sections, we will explore how this affects predictions for 
the Lyman-a forest. 



3.2 Effective optical depths and photoionisation 
rates 

To obtain detailed predictions for the Lyman-a forest in 
our simulations, we create synthetic spectra, as described 
in Sect. 12.31 Examples of such spectra are shown in Fig. [1] 
They are based on neutral hydrogen fractions computed un- 
der the assumption of ionisation equilibrium with the UV 
background which is described by the FG'09 model. In com- 
parison to the simulation without blazar heating, there is 
significantly more flux transmitted in the runs that include 
blazars. Individual absorption lines are shallower and max- 
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Figure 5. We compare IGM temperatures in simulations with 
and without blaza r heating to constraints from Lyman-o; for- 
est observations by I Becker et al.l (|201ll ). The temperatures in the 
simulations were co mputed at th e same redshift-dependent over- 
density A(z) as in I Becker et al.l (|201l|), so as to allow a direct 
comparison. At 2; < 3, the additional heat input in the simula- 
tion with intermediate blazar-heating brings the IGM tempera- 
ture into excellent agreement with the observations. At z > 3, 
where the temperature increase due to blazar heating becomes 
progressively subdominant, the IGM is too hot in all our simula- 
tions, including the one without blazar heating. This discrepancy 
indicates tension between the He ll reionisation history implied by 
the FG'09 UV background and the Lyman-ck; forest temperature 
measurements. 
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Figure 6. Effective optical depths, Tgff , as a function of redshift 
for simulations w ith and wi t hout blazar heating. Ob s ervati onal 
data from FG'08. IViel et al.l (l2004l V and lTvtler et al.l (l2004l ) are 
shown for reference. The simulations with blazar heating repro- 
duce the observed optical depths very well. The results are based 
on the self- consistent UV background. 



imum values close to complete transmission are attained. 
Blazar heating boosts the transmitted flux most efficiently 
in lines with intermediate column densities, while fully ab- 
sorbed lines are almost unaffected. Evidently, due to these 
substantial differences, the mean transmission in one of 
the simulations will be inconsistent with observational con- 
straints. 

In the past, simulated absorption spectra for the 
Lyman-a forest were typically not taken at face value. In- 
stead, the optical depths were rescaled to match the ob- 
served mean transmission values at every redshift separately. 
This can be justified if the actual UV background in the 
Universe is only poorly constrained. We checked that chang- 
ing the UV background does not have a noticeable dynami- 
cal impact in the simulations. However, it will significantly 
change neutral hydrogen fractions, which are computed un- 
der the assumption of ionisation equilibrium in the presence 
of this background. In short, increasing the UV flux by some 
factor boosts the photoionisation rate by the same factor, 
and thereby lowers the lifetime of a freshly recombined hy- 
drogen atom also by that factor. In equilibrium, the neutral 
hydrogen density (and thus also the optical depths) are pro- 
portional to the lifetime of a neutral hydrogen atom. Thus, 
by tuning the UV background one can rescale the optical 
depths. Commonly this was applied in the following way: 
using densities and internal energies from simulations, the 
corresponding optical depths were computed. These were 
then rescaled by an appropriate change of the UV back- 
ground such that the resulting mean transmission was in 
agreement with observations. 

In the following, we check whether such a rescaling of 
the UV background is necessary to make our simulations 
match observed mean transmission values. Figure [6] displays 
the effective optical depth, Teff, corresponding to the mean 
transmission defined by e"^^^^ = {^~^)^ ^s a function of red- 
shift. Our simulations with blazar heating nicely reproduce 
the observed effective optical depths at all considered red- 
shifts, while the simulation without blazar heating yields val- 
ues that are significantly too large. This shows that blazar 
heating together with the recent FG'09 model of the cos- 
mic UV background can fully account for the observed mean 
transmissions as a function of redshift. Instead, the run with- 
out blazar heating requires a rescaling of the UV background 
to get mean transmissions that are compatible with the data. 

We perform our analysis of the simulated Lyman-a for- 
est in the remainder of this paper in two independent ways, 
these are: 

• We compute synthetic Lyman-a forest absorption spec- 
tra by self-consistently adopting the same FG'09 UV back- 
ground model that was used to follow photoheating in the 
simulations. This UV background model agrees well with ob- 
servational constraints from FG'08. In the following we refer 
to results obtained in this way to being based on the self- 
consistent UV background In this approach, all simulations 
are treated equally and the change in the Lyman-a forest 



^ Self- consistent here refers to the fact that the same UV back- 
ground model is used for the simulations and the Lyman-a forest 
analysis. The model is not based on following the UV sources in 
the simulations. 
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properties due to blazar heating is not absorbed by artifi- 
cally rescaling the mean transmission to the observed value. 
We win later show that this UV background allows us to si- 
multaneously match observed IGM temperatures, Lyman-a 
forest line widths, and mean transmission values. 

• In an alternative analysis, we artificially rescale the UV 
background at every redshift such that the optical depths of 
our simulated absorption spectra match the observed mean 
transmission. We, thus, compare the different simulations 
at the same mean transmission value. Results obtained in 
this way will be referred to being based on the tuned UV 
background. 

Overall, these results show that in the run without 
blazar heating realisitic optical depths can only be obtained 
by rescaling the UV background. In the runs with blazar 
heating, the measured mean transmissions is, instead, al- 
ready almost perfectly matched for all redshifts 1.9 < z < 3 
when using the UV background that was adopted in the sim- 
ulations. Thus, results obtained by the self- consistent and 
the tuned analysis differ only very mildly for these runs. 



3.3 Transmitted flux probability distribution 
functions and power spectra 

We now move on to consider the PDFs and power spec- 
tra of the transmitted flux fraction at different redshifts. 
Both statistics are s ensitive to the thermal state of the IGM 
([Theuns et al.ll2000l ) . Before analysing our synthetic spectra 
we smooth them to a resolution that is typical for real obser- 
vations and add appropriate noise that accounts for Poisson 
distributed photon shot noise as well as a constant readout 
noise. We choose a Gaussian filter with a = 3.5kms~^ and 
assume a signal-to-noise ratio that corresponds to a value of 
60 in 2.5kms~^ pixels. In addition, readout noise is added 
with a standard deviation corresponding to 0.5% of the con- 
tinuum flux. 

Figure [7] displays transmitted flux PDFs for simulations 
with and witho ut blazar heating and compares them to ob- 
servations from lKim et al.l (120071 ). For our self- consistent UV 
background, we find that in the intermediate transmitted 
flux range F = 0.05 — 0.9, blazar heating substantially lowers 
the PDF, so that it matches the data. Close to full transmis- 
sion the behaviour is different. At z = 2.94, the PDF from 
the simulation without blazars drops near F = 1, which is 
not observed in the data. However, also the simulations with 
blazar heating deviate somewhat from the data at F = 1. 
One should keep in mind though that the observed PDFs 
are very sensitive to uncertainties in the adopted contin- 
uum level there. An error of a few percent tends to broaden 
and lower sharp peaks in the PDF near F = 1, like those 
found in the runs with blazar heating. This might account 
for a large part of the residual deviation between the data 
and our blazar runs. A discussion of the unc e rtaint ies in the 
PDF near F = 1 can be found in lViel et aJ] (|2009l ). 

We now compare flux PDFs after rescaling all sim- 
ulations to the observe d mean tr a nsmis sion of the corre- 
sponding data sets by iKim et al.l (|2007l ), i.e. we use the 
tuned UV background. As expected, we find that the differ- 
ences between the different simulations are reduced. How- 
ever, the runs with blazar heating are still in significantly 
better agreement with the observations. In particular, they 



do not overproduce the flux PDF in the F — 0.5 — 0.9 trans- 
mitted flux range. iBolton et al.l (|2008l l showed that lower 
values of the flux PDF in that range result from a flat- 
ter or inverted T-p relation, which is naturally produced by 
blazar heating. In contrast, higher IGM temperatures alone 
are not sufficient to match observations there. Especially at 
z — 2.07, part of the discrepancy in the run without blazar 
heating may, how ever, be due to finite box size effects (see 
iTvtler et al.ll2009l V At higher redshift, where the deviation is 
larger and finite box size effects are expected to be smaller, 
this is in contrast unlikely to significantly affect the compar- 
ison. 

As a next step, we compute one-dimensional 
power spectra of the transmitted flux contrast, 
exp(— r)/(exp(— r)) — 1, along our lines of sight. They 
are shown for runs with and without blazar heating at two 
different redshifts in Figu re [51 and are com pare d to obser- 
vationa l constraints from IViel et al.l ([200J) and iKim et al.l 
(|2004al V Based on the self- consistent UV background, we 
find that blazar heating drastically improves the agreement 
between simulated and measured values, reproducing both 
the shape and the normalisation of the observed power 
spectra nicely. The only remaining small deviation is in the 
z — 2.52 simulation results at large k values. At these small 
scales the measured power spectrum is, however, sensitive 
to residual contaminatio n by metal lines that can result in 
excess power there fe.g.. iMcDonald et al.ll2000l V 

When comparing power spectra based on our tuned UV 
background analysis, we rescale optical depths to match the 
ob served mea n trans missi ons of the correspo nding data sets 
bv I Viel et all (|2004h and iKim et aD (l2004al ). We find that 
the runs with blazar heating are in significantly better agree- 
ment with the observations on large scalesfj On intermedi- 
ate and small scales, where the difference caused by blazar 
heating is smaller for rescaled optical depths, all simulation 
models provide an acceptable match to the data. 

To summarise, both the one and two-point statistics of 
the transmitted flux are in excellent agreement with obser- 
vations once blazar heating is accounted for in simulation 
predictions, independent of whether or not we allow for a 
rescaling of the UV background to match the observed mean 
transmissions. 



3.4 Voigt proflle fltting of Lyman-o; lines 

Another technique that has often been used to analyse ob- 
served Lyman-a forest spectra is to decompose them into in- 
dividual absorption line s by fitting theni with superpositions 
of V oigt profiles (e.g., IHu et all [19951 : iKirkman k Tvtlerl 
[19971) . Voigt profiles describe a thermally broadened spectral 
line of finite lifetime, in other words, they can be obtained 
by convolving a Lorentzian line profile with a Gaussian. The 
thermal broadening is usually quantified using the broad- 
ening parameter b = (2/cT/mH)"^^^, where k is the Boltz- 
mann constant, T the temperature, and mn the mass of 
an hydrogen atom. Exploiting this relation, the distribu- 
tion of line widths can be linked to the thermal state of the 



^ Finite box-size effects are not expected to singificantly affect 
our simulated power spectra on sca les k > 0.005 skm"-*^ (see, e.g., 
iTvtler et al]|2009l : IViel et al.|[20ol) . 
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Figure 7. PDFs of the transmitted flux fraction for simulations with and without blazar heating. Results are shown for redshifts z = 2.07 , 
2.52, and 2.94, and for the three different normalisations of the blazar heating rate. Observational constraints from iKim et al.l ([20071) 
are shown for reference. The left panels show results for the self- consistent FG'09 UV background, while the right panels show results 
for a UV background that, for each simulation and redshift, was tuned to match the observed mean transmission. In both cases, the 
simulations with blazar heating are in substantially better agreement with the data than the simulations without it. 



IGM. Thus, measuring the former at diflPeren t redshifts al- 
lows tracing the therma l history of the IGM (|Schave et al.l 
I2OOOI : iRicotti et alllioool ). 

We applied this method to our synthetic Lyman-a for- 
est spectra and compared the distributions of line- widths, 
b, and neutral hydrogen column densities, A^hi, to observa- 
tions. The Voigt-p rofile-decompositi on was performed with 
the AuTOVP code (JDave et al.lll997l V Before applying it to 
our synthetic spectra we smoothed them and added noise as 
described in Sect. [331 

Figure [9] shows the distribution of absorption lines in 
b and A^hi for simulations without and with intermediate 
blazar hea ting. The lower 6 -cutoff , obtained from observa- 
tions by Ki rkman fc Tvtlerl (| 19971 ) , is shown for reference as 
a function of column density. As expected, the blazar heating 



boosts the IGM temperatures, which results in larger ther- 
mal broadening and accordingly shifts the lower 6-cutoflP to 
higher values. This brings the simulation prediction into ex- 
cellent agreement with the measured cutoflF if blazar heating 
is included. On the other hand, the run without blazar heat- 
ing yields absorption lines with 6- values that are too small. 
A direct comparison of the b-Nui distribution of the ob- 
served sample and synthetic samples of similar size is given 
in Appendix [A] This comparison corroborates the superb 
agreement between the data and our intermediate blazar 
heating model. 

The shift of the lower ^-cutoflF can also be directly seen 
in the normalised distribution function of 6- values, which 
is given in Figure [10] for lines w ith A^hi > lO^^cm-^. Ob - 
servational constraints based on iKirkman fc Tvtlerl (| 19971 ) 
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Figure 8. One-dimensional power spectra of transmitted flux contrast exp(— r)/(exp(— r)} — 1. Results are sh own f or simulations with 
and without blazar heating for redshifts z = 2.07 and 2.52. Observational constraints from I Viel et al.l (|200J) and iKim et al.l (|2004al ) 
for redshifts z = 2.125 and 2.58, respectively, are shown for reference. The left panels show results for the self- consistent FG'09 UV 
background, while the right panels show results for a UV background that, for each simulation and redshift, was tuned to match the 
observed mean transmission. Including blazar heating in the simulations generally improves the agreement with this data, in particular 
for our self- consistent UV background analysis. 



are shown for reference. For both the self- consistent and the 
tuned UV background analysis, the data are only consistent 
with simulations that include blazar heating. Comparing the 
figure's panels also shows that this method relies only very 
weakly on the assumed photoionisation rate, as a rescaling 
of the optical depths only changes the column densities but 
not the widths of absorption lines. Thus, the only effect is 
a slight shift of the threshold value above which lines are 
included. Due to the weak dependence of the lower 5-cutoff 
on column density, this does not strongly affect the results. 

Figure [11] compares both, the normalised distribution 
of column densities, as well as the frequency of lines with 
given column density per unit redshift to observations. The 
former distribution is only mildly affected by blazar heat- 
ingfl, and all simulations are in good agreement with the 
data. In contrast, the frequency of lines per unit redshift 
is much better reproduced in the runs with blazar heating. 
This is, of course, closely related to the better agreement of 
flux PDFs in these runs. 

In sum, we find that the simulation predictions match 
the observed properties of Lyman-a absorption systems im- 



^ The distributions depend however signiflcantly on the choice of 
the lower column density cutoff. 



pressively well when blazar heating is included. In partic- 
ular, we would like to stress that the improved agreement 
with the observed lower cutoff of line- widths is a direct and 
independent confirmation of high IGM temperatures. 



4 OUTLOOK AND COSMOLOGICAL 
IMPLICATIONS 

In this section, we use our previous results on the effects 
of blazar heating and qualitatively discuss their impact 
on the local Lyman-a forest, on estimates of ag (through 
the Lyman-a forest), on constraints of (sterile) neutrino 
masses, Helium ii reionisation, formation of dwarfs and 
galaxy groups/clusters, and on particle acceleration at cos- 
mological structure formation shocks. While we provide per- 
tinent arguments for each point, we stress that each of these 
topics merits a thorough analyses to reliably quantify and 
confirm the mentioned effects, something which is beyond 
the scope of the present work. 

The trend that the Lyman-a forest probes larger over- 
densities at lower redshift, as mentioned in the discussio n of 
Fig. O continues down to z = (see e.g., ISchavdbOOlh . In 
particular, significant absorption lines in the local Lyman-a 
forest are expected to be due to absorbers at overdensities 
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Figure 11. The normalised distribution of column densities A^hi of Lyman-o; lines in simulations with and without blazar heatin g at 
redshift z = 3 {top panels). Included are all lines with A^hi > lO-'^^cm"^. Observational constraints based on Kir kman &l Tvtlerl (|l997) are 
given for reference. They were obtained by using all lines from their Table 1 which were classified as Lyman-a lines, have A^hi > lO-'^^cm"^ 
and redshifts close to the simulation output (2.75 < z < 3.05). The lo wer panels shows the number of lines per unit redshift and per dex 
in column density. Observational data from lKirkman fc Tvtlerl (|l997[ ) are shown for comparison, where only lines with 2.88 < z < 3.05 
were used to exclude gaps in their spectra at other redshifts. Results are shown for the self- consistent UV background (left panels) and 
the tuned UV background (right panels). All error bars assume Poissonian errors. 



5 ^ A ^ 100 (|Dave et al.ll20ld ). where the impact of blazar 
heating is strongly reduced. For example at an overdensity 
of 10 and 2; = 0, we find only a ^ 30% temperature differ- 
ence due to blazar heating instead of the factor of 10 dif- 
ference that is shown in Fig. [2] for mean density. We, thus, 
do not expect very drastic effects of blazar heating on the 
local Lyman-a forest. We would also like to point out that 
the latter does not yield direct constraints on the thermal 
state of the IGM at mean density. In particular, power-law 
extrapolations of the measured T-p relation down to the 
mean density fe.g lRicotti et al.ll2000l ) are invalid in a blazar- 
heated universe in which the T-p relation is tilted. This is 
illustrated in Fig. 1121 It may be interesting to explore in 
detail how blazar heating changes the statistics of the local 
Lyman-a forest in a future work. 

One interesting aspect of our results is that by including 



blazar heating we were able to reproduce the observed trans- 
mitted flux power spectra well with simulations that adopt 
as = 0.809. This suggests that blazar heating could alleviate 
reported tensions between as-values inferred from Lyman-a 
forest studies and other cosm olo^ical constraints, mo st no- 
tably CMB observations fe.g.. lViel fc HaehneltllJOOel ). This 
tension is also present in the flux power spectra in our tuned 
UV background analysis (see right panels of Fig. [8]). There, 
the large-scale flux power in simulations without blazar heat- 
ing falls short of the observed level. A joined analysis of 
flux PDFs of 27 high-resolution, high signal-to-noise ratio 
Lyman-a forest spectra (the Large UVES Quasar Ab sorp- 
tion Spectra, or LUQAS, sample bv lKim et al.ll2004br ) and 
the SDSS flux power spectrum which is based on a large 
number of quasar s pectra (3035) with low -resolution and 
low signal- to- noise ([McDonald et al.l l2005l l yields a larger 



© 0000 RAS, MNRAS 000, 000-000 



14 E. Puchwein, C. Pfrommer, V. Springel, A. E. Broderick & P. Chang 



no blazar heating 



intermediate blazar heating 




3 -2 

logio(p/{P» 

Figure 12. The gas density-temperature distribution in simulations without blazar heating (left panel) and with intermediate blazar 
heating (right panel) at redshift z = 0. The gas mass per pixel Mp i y is c olour-coded. The colour scale spans 3 orders of magnitude. 
Also shown is the power-law T-p relation obtained by iRicotti et al.l (|2000) from observations of the local Lyman-o; forest that probes 
overdensities 5 ^ A ^ 100 (shown by the solid line). The dashed line indicates the extrapolation to the mean density (grey circle) 
assuming the same power-law relation holds there. The green diamond indicates the median gas temperatures at mean density, which 
is the same quantity that is shown in Fig. [2] as a function of redshift. The thermal state of the IGM in the range where it is probed by 
the local Lyman-o; forest is fairly similar in runs with and without blazar heating. However, while the power-law extrapolation to mean 
density works remarkably well without blazar heating, it completely fails in a blazar-heated universe in which the T-p relation is tilted 
upwards at low density. 



value of as = 0.86 di 0.03 (|Viel et al.ll2009l ) than an analysis 
of CMB data alone. This discrepancy appears to be related 
to the higher temperature To that the flux power spectrum 
prefers and the fiducial T-p relation slope of 7 ^ 1 around 
which the model was expanded. In fact, there is a degen- 
eracy between the logarithmic slope of the T-p relation at 
mean density, 7 = dlogP/dlogp, and as; larger values of 7 
imply increasing values for erg to be consistent with the data 
(see Figure 2 in lViel et al.ll2009l V For a T-p relation with a 
positive slope, 7 > 1, and fixing all other parameters, the 
inferred value of as > 0.9 would be in strong disagreement 
with values inferred from the CMB or other cosmological 
probes. (We note, however, that such a procedure of fixing 
all but two parameters is only for illustrative purposes; in 
practice, all necessary parameters have to be allowed to vary 
simultaneously in order to identify the global best fit to the 
data.) Blazar heating provides a physically motivated sce- 
nario for an inverted T-p relation with 7 2^ 0.5, allowing 
the normalisation to match current concordance values for 
as — 0.8 within la error bars. 

Previous works that analysed Lyman-a power spectra 
in order to derive cosmological parameters and constraints 
on the mass of a sterile neutrino imposed Bayesian priors o n 
the T-p relation (JMcDonald et al.ll2005l : ISeliak et al.ll2006l ). 
Those were primarily m otivated by theoretic ally expected 
limits, < 7 - 1 < 0.6 (JHui fc Gnedinlll997l ). rather than 
informed by Lyman-a forest observations themselves. These 



adopted priors in the parameter estimation algorithm ex- 
cluded the possibility of an inverted T-p relation during the 
marginalisation process. Hence, this may have lead to biased 
results if this shortcoming has not been absorbed somehow 
by the other (empirical) parameters within the minimisa- 
tion procedure (which employed a total of 34 parameters in 
iMcDonald et al.ll2005l ). We have shown in Figure [8] that the 
cutoff scale of the one-dimensional power spectrum of the 
transmitted flux contrast is severely modified as a result of 
blazar heating and shifted towards lower values in compar- 
ison to models with photoheating alone. In principle, this 
could allow for a larger intrinsic cutoff in the transfer func- 
tion as potentially imposed by a warm dark matter model 
and, hence, might weaken the constraints on the neutrino 
masses. 

So far. He 11 reionisation seemed to be the only pro- 
cess able to substantially contribute to the late-time heat- 
ing of the IGM. As a result, inhomogeneous He 11 reionisa- 
tion has been studied to explain the high temperatures and 
the inverted T-p relation of the IGM inferred fr om Lyman- 
a forest measurements (|Furlanetto Sz Ohll2008l ). However, 
some of the properties that have been attributed to He 11 
reionisation need to be re-evaluated within the framework of 
blazar heating, which can account for substantial heating for 
z < 3.5. It will be an interesting challenge to constrain both, 
the characteristics of He 11 reionisation and blazar heating. 
If successful, this might even allow a determination of the 
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redshift evolution of the blazar luminosity densit y by em- 
ploying an analogous argument to that by Soltan| (|l982l l for 
the thermal history of the IGM. This could be achieved by 
combining the detailed temperature history as inferred from 
statistical analyses of the Lyman-a forest fe.g. jBecker et alj 
1 20 111 ) with the quasar luminosity density and improved con- 
straints on the intrinsic quasar spectra above 4 Ry, the en- 
ergy range that is of importance for He ii reionisation, as a 
function of redshift. 

A dramatically increasing minimum entropy of the IGM 
for z < 3 has a number of important impli cations for struc - 
ture formation (discussed in detail in iPfrommer et alfcOllh . 
Such a redshift dependent entropy floor increases the char- 
acteristic halo mass below which dwarf galaxies cannot form 
by a factor of ^ 10 at mean density over that found in mod- 
els that include photoionisation alone (when adopting our 
successful intermediate blazar heating model which is only 
so mewhat less efficient than the optimistic model studied 
in lPfrommer et al.ll201lh . This may help resolve the "miss- 
ing satellite problem" in the Milky Way and the "void phe- 
nomenon" of the low observed abundances of dwarf satellites 
compared to cold dark matter simulations. 

Simultaneously, the entropy injection by blazars is in 
some sense an amalgam of both the preheating and AGN 
feedback mechanisms for galaxy clusters. However, blazar 
heating operates on much larger distances than AGN feed- 
back and provides a tim^e- dependent entropy injection rate, 
peaking near 2; '-^ 1 (in contrast to the instantaneous heating 
previously assumed for preheating models). This suggests a 
scenario for the origin of the cool core/non-cool core bi- 
modality in galaxy clusters and groups where early forming 
galaxy groups {z > 1) are unaffected due to the small level 
of entropy in the IGM at these times that can be quickly 
radiated away. However, late forming groups do not have 
sufficient time to cool before the entropy is gravitationally 
reprocessed through successive mergers. 

In Figure 2] we have shown the redshift evolution of 
the volume-weighted PDF of density in the blazar heat- 
ing model, which starts to deviate from the model without 
blazar heating already ai z ^ 2. By 2; < 1, the PDF is no- 
ticeably narrower in the blazar heating model, and the wings 
at low and high densities are suppressed. This suggests that 
blazar heating is indeed able to slow down non-linear struc- 
ture formation on these scales. The degree to which this 
affects the luminosity function of dwarfs and the bimodal- 
ity of central entropy values in galaxy clusters is subject of 
hydrodynamical cosmological simulations that we plan to 
carry out in forthcoming work. 

Additionally, we have show in Fig. 2] that the median 
of the volume- weighted temperature PDF increases by more 
than a factor of 100 at 2; = as a result of blazar heating. 
This implies an increase of the sound speed, Cs oc vkT, by 
a factor of ^ 10. Hence the strengths of the first shocks (as 
measured by the Mach number, M = t'shock/cs) which track 
the collapse of the IGM from voids onto large scale struc- 
ture filaments an d super-cluster regions are also reduced by 
a factor of ^ 10 (|Pfrommer et al.ll2006l V This should have 
important implications for diffusive particle acceleration at 
structure formation shocks as the acceler ation efficiency is a 
strong function of M ach number (see e.g. jEnBlin et ahlbOOTI : 
iKang Sz Rvull2010l , and references therein). Similarly, this 
could also have interesting consequences for models of the 



resulting non-thermal emission from galaxy clusters, such as 
diffuse radio halos and relics. 



5 SUMMARY AND CONCLUSIONS 

Recent studies discovered that TeV 7-ray emission from 
blazars can significantly heat the low-density IGM by anni- 
hilating and pair producing on the extragalactic background 
light. The resulting high-energy electrons and positrons ex- 
cite plasma be am instabilities which dissipate the pairs ' 
energy locally (iBroderick et al.1 I2OIII : IChang et aP l201lh . 
leading to important implications for structure formation 
(|Pfrommer et al.ll201lh . In this picture, the IGM acts as a 
cosmic heat bath that absorbs most of the bolometric energy 
of TeV blazars, which are ultimately powered by accretion 
onto their central super-massive black holes. Since the dis- 
tributions of EBL photons and of TeV blazars are nearly 
homogeneous on cosmological scales, the resulting pair den- 
sity is spatially homogeneous as well. As a result, the implied 
heating rate is independent of the IGM density. 

We have included the expected heating term in state- 
of-the-art cosmological hydrodynamical simulations and in- 
vestigated its impact on the thermal history of the IGM and 
the properties of the Lyman-a forest. We find that: 

• In agreement with linear theory predictions by 
IChang et al.l (|201lh , the temperature of the low density IGM 
is strongly boosted by blazar heating, resulting in an in- 
verted temperature-d ensity relation as favoured by Lyman- 
a forest observations (|Viel et al.ll2009l V 

• At z < 3, when blazar heating is the dominant heating 
source, the IGM temperature, as a function of redshift, is 
also in good agreement with current constraints based on the 
curvature of observed Lyman-a forest spectra ([Becker et al.l 

ionl). 

• The mean transmissions predicted by runs with blazar 
heating are in excellent agreement with observations when 
employing the photoionisation rates predicted by up-to-date 
models of the cosmic UV background (FG'09). We adopt 
these rates in our self- consistent analysis which uses the 
same UV background for simulations and analysis. To allow 
for the uncertainties in the determination of the photoion- 
isation rates, we additionally explore how results change if 
we tune the UV background to match the observed mean 
transmission values. 

• Blazar heating alters the distribution functions and 
power spectra of the transmitted flux and resolves discrepan- 
cies with observational constraints that exist when photoion- 
isation is the only modelled source of heating. This improved 
match of simulations and data in our blazar heating models 
is also true if we rescale the photoionisation background to 
match the observed mean transmission at each redshift. This 
seems to be a uniq ue property of the inverted T-p relation 
([Bolton et al.ll2008l ) which is naturally produced by blazar 
heating. 

• Heating of the IGM by blazars shifts the lower cutoff 
of Lyman-a absorption line widths to larger values due to 
increased thermal broadening and decreases the abundance 
of weak lines with column densities A^hi ^ 3 x lO^^cm"^. 
This makes the lower cutoff and the abundance of lines per 
unit redshift and column densities consistent with observa- 
tions and provides additional evidence for high IGM tem- 
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peratures independent of whether or not we rescale to the 
observed mean transmission. 

• The inclusion of blazar heating brings the simulated 
large-scale flux power into agreement with the observed 
level. The simulation predictions are based on the most re- 
cent cosmological parameters which suggests that previous 
tensions in Lyman-a measurements of ag when compared 
to other cosmological constraints are alleviated by this novel 
heating mechanism. Our results also support previous claims 
that the empirical assumption of an inverted temperature- 
density relation for th e IGM a lready relieves the tension in 
the constraints on ag (|Viel et a l. 2009). 

• Blazar heating dramatically alters the volume- weighted 
temperature PDF which implies that the strengths of first 
structure formation shocks during non-linear collapse are 
weakened. Future work is needed to work out the impli- 
cations for the particle acceleration efficiencies and the re- 
sulting non-thermal emission from galaxy clusters. Modifi- 
cations in the density PDF furthermore suggest that blazar 
heating may have interesting effects on structure formation, 
particularly on the scale of dwarf galaxies. 

To summarise, including blazar heating in cosmologi- 
cal simulations results in properties of the IGM and the 
Lyman-a forest that are in excellent agreement with ob- 
servations. In particular, we demonstrate in this work that 
our intermediate blazar heating model in combination with 
the FG'09 UV background represents the first successful, 
physically-motivated model that simultaneously reproduces 
all Lyman-a forest data (flux PDFs, power spectra, tem- 
perature evolution, line width and column density distri- 
butions) without the need to rescale the UV background 
at each redshift. This strengthens the case that this novel 
heating mechanism is indeed operating in our Universe, and 
at the same time it resolves long standing difficulties in our 
understanding of the Lyman-a forest. 

This impressive improvement in modelling the Lyman- 
a forest appears to be a direct consequence of the peculiar 
properties of blazar heating. Its heating rate, which is in- 
dependent of IGM density, naturally produces the inverted 
temp erature-den sity relation that Lyman-a forest data sug- 
gest (|Viel et al.l 12009). The recent and continuous nature 
of the heating generates the redshift evolutions of the IGM 
temperature, the effective optical depth, the flux PDF, and 
flux power spectrum that accurately match observed trends 
for the redshift interval 2 < z < 3 analysed here. Finally, the 
magnitude of the heating rate required by Lyman-a forest 
data approximately coincides with the total energy output 
of TeV blazars (or equivalently ^ 0.2% of that of quasars). 
While all these properties are quite unusual in the context 
of H and He i/ii photoheating and traditional feedback pro- 
cesses, e.g., that of AGN or supernova driven winds, they 
are natural in the TeV-blazar heating mechanism. There is 
hence ample motivation for further studies of this scenario, 
which are also needed to establish in full detail how blazar 
heating affects cosmological constraints and structure for- 
mation over a range of scales. 
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APPENDIX A: COMPARING THE B 
DISTRIBUTION TO OBSERVATIONS 



Nu 



In the following, we directly compare the simulated Lyman- 
a line widths, 6, and neut ral hydrogen co lumn densities, 
A^Hi, to observations by K irkman fc Tvtlerl (|l997). Rather 
than showing the distribution function as in Fig. [9l we 
present scatter plots of b and A^hi values. All panels con- 
tain the same number of significant A^hi > 10^^"^cm~^ lines 
as the observed sample. This facilitates a direct comparison. 
Note, however, that these results are more strongly affected 
by small-number statistics than those in Fig. (9] The figure 
displays results for simulations without blazar heating, re- 
sults for simulations with intermediate blazar heating, as 
well as the observed sample. 

Clearly, the lower 6-cutoff is in much better agreement 
with the observations if blazar heating is included. Interest- 
ingly, the region defined by 12.6 < log^Q(A^Hi cm^) < 13.1 
and 15kms~^ < b < 23kms~^, in which quite few lines 
fall in runs with blazar heating (due to the better statistics 
this can be more easily seen in the right panel of Fig. [9|), is 
also very sparsely populated by the observed sample. This 
suggests that by including blazar heating also more subtle 
features of the b — Nui distribution are better reproduced. 
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